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The thin-shell instability has been named as one process, which can generate entangled structures 
in astrophysical plasma on collisional (fluid) scales. It is driven by a spatially varying imbalance 
between the ram pressure of the inflowing upstream plasma and the downstream’s thermal pressure 
at a non-planar shock. Here we show by means of a particle-in-cell (PIC) simulation that an analogue 
process can destabilize a thin shell formed by two interpenetrating, unmagnetized and collisionless 
plasma clouds. The amplitude of the shell’s spatial modulation grows and saturates after about ten 
inverse proton plasma frequencies, when the shell consists of connected piecewise linear patches. 

PACS numbers: 52.35.Tc 52.65.Rr 52.38.Kd 


A shock wave forms if an obstacle moves through a 
medium at a speed, which exceeds the speed at which 
the medium can convey information. Shocks are ubiqui¬ 
tous in fluids and plasmas. The bow shock that separates 
the solar wind from the Earth’s magnetopause is proba¬ 
bly the best-understood plasma shock mm- Even larger 
shocks can form where energetic plasma outflows collide 
with the interstellar medium. The solar wind termina¬ 
tion shock separates the heliosphere from the interstel¬ 
lar medium [3]. The existence of plasma shocks outside 
of the solar system, such as supernova remnant shocks 
[4], can be inferred from the electromagnetic radiation 
they emanate. A comprehensive review of these non- 
relativistic magnetized shocks is given by Ref. [5]. 

The properties of shocks are determined by those of 
the underlying medium; the most important one being 
the particle’s mean free path. If the mean-free path of 
the particles is small compared to the shock scale, then 
we speak of a fluid- or hydrodynamic shock. 

A fluid shock is a discontinuity that separates a pre¬ 
shock (upstream) fluid from a post-shock (downstream) 
fluid. The velocity component of the upstream flow along 
the shock normal and measured in the reference frame 
of the shock is supersonic. This velocity component is 
reduced to a subsonic value by the shock crossing and the 
released directed flow energy is transformed into internal 
energy. The velocity component of the fluid, which is 
orthogonal to the shock normal, is left unchanged by the 
shock crossing. If the shock is planar and the lateral 
velocity constant, then we can remove the lateral velocity 
by choosing a comoving reference frame. 

Let us assume that the upstream velocity is spatially 
uniform and that the shock front has a sinusoidal mod¬ 
ulation. The direction of the shock normal and, hence, 
the components of the upstream velocity vector along this 
normal and orthogonal to it vary along the shock bound¬ 
ary. We can not find a reference frame, in which the 
lateral velocity vanishes everywhere. The normal compo¬ 
nent of the fluid velocity is reduced by the shock crossing 


while the lateral velocity remains unchanged, which im¬ 
plies a rotation of the fluid velocity vector. The rotation 
angle varies with the position along the shock bound¬ 
ary. This spatially varying rotation of the velocity vec¬ 
tor implies a spatially varying momentum transfer from 
the upstream fluid to the shock boundary and, hence, a 
spatially varying imbalance between the ram pressure of 
the upstream flow and the thermal pressure of the down¬ 
stream fluid. Such an imbalance can yield a linear [6] or 
a nonlinear [7] instability in the case of a thin fluid shell 
that is enclosed by two shocks. The thin-shell instability 
(TSI) amplifies the corrugation of the shock boundaries. 
The nonlinear TSI discussed in Ref. [7] has been exam¬ 
ined by means of hydrodynamic [8] and magnetohydro¬ 
dynamic [9] simulations. 

Hydrodynamic pm [TT] and particle-in-cell (PIC) [12]- 
[T4] simulations demonstrate that the ablation of a target 
by an intense laser pulse generates an energetic plasma 
flow. The collision of this plasma flow with ionized resid¬ 
ual gas, which is usually collisionless, gives rise to the 
electrostatic shocks that have been detected in many sim¬ 
ulations and experiments [151123] . The laser-plasma ex¬ 
periment performed in Ref. [23] resulted in a thin plasma 
shell, which was enclosed by two plasma shocks. A laser- 
plasma experiment can thus reproduce a configuration 
that is similar to that discussed in Refs. [7], which gives 
rise to the nonlinear TSI. 

Collisionless shocks are mediated by collective electro¬ 
magnetic fields rather than by binary collisions between 
particles. The nature of the electromagnetic field dis¬ 
tribution and the structure of the shock transition layer 
strongly depend on plasma parameters such as the mag¬ 
netization and the ratio between the shock speed and 
the light speed. A nonrelativistic plasma shock in an un¬ 
magnetized collisionless plasma is mediated by an elec¬ 
trostatic field [TanH]. The source of this electric field is 
the thermal diffusion of electrons from the dense down¬ 
stream plasma into the upstream plasma, which leaves 
behind a positive net charge in the downstream region 
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FIG. 1: (Color online) (a) The collision of the equally dense 
plasma 1 (at rest) and plasma 2 (speed r'o) yields two plasma 
shocks (blue sinusoidally varying curves), which enclose a 
plasma shell. Each shock is displaced along x with respect to 
its average plane (horizontal dashed black lines). The electric 
field E, which mediates the plasma shock, points along the 
shock normal. Panel (b) focuses on the center of (a), where 
the horizontal electric field component Ey is strongest. Ions 
(p) and electrons (e) are deflected into opposite directions as 
they cross a shock and the resulting net current J yields a 
magnetic B-field. 


close to the shock. An electric field builds up that accel¬ 
erates electrons from the upstream region into the down¬ 
stream region. An equilibrium is established when the 
net outflow of electrons from the downstream region due 
to thermal diffusion is balanced by the electrons that are 
injected back into the downstream region by the elec¬ 
tric field. The electric field is antiparallel to the density 
gradient. 

The electric field also affects the ions. It slows down 
the ions along the shock normal when they cross the 
shock and move into the downstream region. Such a 
structure resembles a shock in a collisional fluid and it 
is thus called an electrostatic shock. In contrast to fluid 
shocks, electrostatic shocks can reflect some of the in¬ 
coming upstream ions and they do not fully thermalize 
the ions that cross the shock. Some of the downstream 
ions are accelerated upstream by the ambipolar electric 
field and such a structure is known as a double layer. A 
plasma shock is usually a combination of a double layer 
and of an electrostatic shock [24l [25] . 

The question arises whether or not an instability, which 
is similar to the TSI, can develop in a collisionless lab¬ 
oratory plasma. Here we study this scenario by means 
of a particle-in-cell (PIC) simulation [26] study with the 
EPOCH code [27]. Consider the case depicted in Fig. 
[^a), where an upstream medium (plasma 2) is moving 
at a uniform flow speed towards the shock (lower blue 
sine curve). The shock boundary is displaced along this 
direction and the amplitude of the displacement varies 
sinusoidally with y. The same holds for the electric field 
vector E that mediates this plasma shock. 

The direction of E varies with y, while the ram pres¬ 


sure force of the upstream flow points along its uniform 
flow velocity vector; the latter is anti-parallel to the av¬ 
erage shock normal. Hence the thermal pressure does 
not cancel out the ram pressure everywhere. The fluids 
are deflected into opposite horizontal directions by both 
shocks along a vertical cut in Fig. Ba) and the deflec¬ 
tion angle increases with an increasing obliquity of the 
shock boundary. Some important consequences of this 
deflection are revealed by our simulation study. 

The plasma shocks are created in the simulation by the 
collision of two plasma clouds, each consisting of protons 
and electrons with the mass ratio mpjme = 1836. The 
electrons and protons of each cloud have an equal density 
no and mean speed modulus. Each cloud is thus initially 
charge- and current-neutral and we can set the electric 
E-field and the magnetic B-field to zero at the time t = 0. 

Plasma 1 is at rest in the simulation frame and the 
mean speed of plasma 2 is i;o. The electron temperature 
Tg = 1.74 X 10^ K and the proton temperature Tp = Tejb 
are representative for laser-produced plasma. Initially, 
the plasma clouds are spatially separated and they are in 
contact at a boundary xsiv) with a sinusoidally varying 
displacement from a straight average boundary. Plasma 
2 moves into plasma 1 for increasing times t > 0 and a 
thin shell forms, which resembles that depicted in Fig. 
Ufa). The ions slow down as they cross the plasma shock 
and their piling up increases the density of the thin shell. 
The density increase depends on the Mach number of 
the cloud velocity with respect to the ion acoustic speed 
Cs = ilck B{TeTp)/mp) ■ {kB- Boltzmann’s constant) 
and measured in the reference frame of the shell. Typical 
values range from just over 2no to 3no for slow planar 
plasma shocks [T71I1H]. The ion acoustic speed is Cs ~ 
5.4 X 10^ m/s for our plasma parameters and the collision 
speed is vq = 3.8cs. We have assumed that the adiabatic 
index is 7 c = 5/3 in order to obtain a reference value for 
Cs. The degrees of freedom in a plasma and hence the 
adiabatic constant depend on the field distribution. 

1 / ^ 

The proton Debye length Xd = {t^^kBTpjn^^^y 
normalizes space. The proton plasma frequency bjp = 

(noe^/eomp) ^ normalizes time. The electric fields are 
normalized to mpVth^pl^ and the magnetic fields to 

mpbjpje^ where Vth = is the protons’ ther¬ 

mal speed. The densities are expressed in units of no- 

The simulation resolves the interval —864 < x < 864 
along the collision direction and 0 < ^ < 180 perpen¬ 
dicular to it. The boundary condition at x = —864 is 
open, that at x = 864 is reflecting and the boundary 
conditions along y are periodic. The collision boundary 
is defined by XB{y) = 1-44 • sin(27r^/90) and it is thus 
centered around x = 0. Plasma 1 initially occupies the 
interval XB{y) < x < 864. Plasma 2 is located in the 
interval —864 < x < XB{y)- Momentum conservation 
implies that the shell will propagate to increasing values 
of X with the speed ^ vq/ 2. The simulation grid resolves 
X by 2400 cells and y by 250 cells. The side length of each 
cell is thus 0.72Xd. Electrons and protons are each rep- 
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resented by 200 particles per cell at the simulation’s start 
and we do not inject new particles while the simulation 
is running. The total simulation time is tgim = 68. 

Figure [^displays the total proton density distribution 
np{x,y) together with Ex{x,y), Ey{x^y) and Bz{x^y) at 
the time t = 4.5. Movie 1 m shows the time-evolution 
of np{x^ y) throughout the simulation and Movie 2 shows 
that of Bz{xB) [32]. The shell is characterized in Fig. 



FIG. 2: (Color online) The thin shell at the time t — 4.5: 
Panel (a) show np{x,y). Panels (b,c,d) show the Ex, Ey and 
Bz components, respectively. The value of Bz has been mul¬ 
tiplied by the factor 10^. 

I^a) by a peak proton density ^2.3, which implies that 
the inflowing ions have been slowed down and compressed 
by the positive potential of the shell. The displacement 
along X of the shell is comparable to the initial amplitude 
of xsiy) at this time. The in-plane electric field vector 
{Ex,Ey) varies with y and the peak amplitude of Ey is 
about one half of that of E^. 

The inflowing electrons and protons are deflected into 
opposite directions when they cross a plasma shock with 
its strong associated F^^-fleld. The therefrom resulting 
net current in the simulation plane generates the ob¬ 
served Bz 0 within the thin shell. We can see this 
from Fig. Bb), which corresponds to a position ^ = 45 
in Fig. The current, which is generated by the hori¬ 
zontal deflection of the electrons and protons, points to 
the left at the lower plasma shock and to the right at the 
upper plasma shock at ^ = 45. Such a current distribu¬ 
tion yields a magnetic field, which points into the plot’s 
plane and, thus, into the same direction as the positive 
z-axis in Fig. Bd). The periodicity of the distributions 
of Ey and Bz equals that of the initial boundary oscil¬ 
lation xsiy) and their phases are shifted by 7r/2, which 
is further evidence for a magnetic field generation by the 
particle deflection by Ey. 

Figure shows the simulation results at tsim = 24. 
The shell has expanded along x and shows a sinusoidal 
displacement with the same phase and period as that 



FIG. 3: (Golor online) The thin shell at the time t = 24: 
Panel (a) show np{x,y). Panels (b,c,d) show the Ex, Ey and 
Bz components, respectively. The value of Bz has been mul¬ 
tiplied by 10^ 


of XB{y)‘ The amplitude of the oscillations of the thin 
shell along x is about 20, which exceeds that of xsiy) 
by a factor 15. The proton density peaks at extrema of 
the oscillations of the thin shell. The strongest struc¬ 
tures in the electric field in Fig. [^b) outline the proton 
density jump between the thin shell and the surrounding 
plasma and they demarkate the locations of both plasma 
shocks. The Ey component in Fig. j^c) shows piecewise 
linear patches that are correlated with the locations of 
the plasma shocks. The peak amplitudes of Ex and Ey 
are comparable in strength. The extent along y of the 
patches in Fig. [^c) is about 30 and their position along 
X varies by about the same value. The angle between 
the electric field {Ex, Ey) and the average normal of each 
shock is thus about 7r/2 within these intervals. The mag¬ 
netic field structures in Fig. |^d) are strongest within 
the thin shell and they remain closely correlated with 
the electric Ey field patches. The numerical values of 
Bz correspond to the ratio between the proton cyclotron 
frequency and the proton plasma frequency and the peak 
value is about 5 x 10“^. The effects of the magnetic field 
on the protons is negligible. 

Figure H compares the electric field contour 
{E^ + Ey)^ = 0.23 with the proton density dis¬ 
tribution of plasma 2. These protons flow towards 
increasing values of x. The proton density at low values 
of X is close to 1. These protons correspond to the 
upstream plasma, before it has encountered the lower 
plasma shock. The lower shock is outlined by the 
lower pair of countour lines. The proton density hardly 
increases when the protons cross the convex parts of 
the lower shock at ^ ~ 65 or ^ ^ 155. The density 
is increased substantially when the protons cross the 
concave parts of the lower shock at y ^ 25 or y ^ 110. 
The density changes gradually as a function of y in 
the intervals where the lower shock has a constant tilt 
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FIG. 4: (Color online) The density distribution of the protons 
of plasma 2 at the time t = 24. Overplotted is the contour 
0.23 of the electric field modulus {E^ . 


angle relative to the incoming upstream plasma flow 
direction. The density change is symmetric with respect 
to the locations of the extrema of the shell’s oscillation. 
The density change is also closely correlated with the 
electric field contours. The proton distribution evidences 
that the inflowing protons are deflected by the electric 
F^y-component towards those y-intervals, where the 
lower shock is concave. The protons accumulate at this 
location and their density increases to a value ~ 1.7. 
This proton flow direction matches that in Fig. [^b). 

We find smaller peaks of the proton density above the 
upper plasma shock, which is outlined by the upper pair 
of electric field contour lines. The protons, which have 
crossed the thin shell and have reached its opposite side, 
are reaccelerated by the upper shock and escape into 
the upper upstream region. This reacceleration is ac¬ 
complished by the double layer component of the upper 
plasma shock. The tilted electric field channels these re¬ 
accelerated protons into the areas centered at y ^ 6b and 
X ~ 120 as well as ^ ~ 155 and x ~ 120. 

The steady deflection of inflowing protons by the hor¬ 
izontal electric field component towards the concave 
shocks will increase further their density at these loca¬ 
tions; the proton distribution in the Figs, [^andj^is not 
in a steady state (See also movie 1). However, the proton 
density difference between the position {x^y) = (20,115) 
and (x^y) = (45,95) in Fig. [^a) is already about 0.5 
and, thus, comparable to the difference between the up¬ 
stream plasma and the downstream plasma. The am- 
bipolar electric field, which is associated with the proton 
density variations along the thin shell, will eventually be¬ 
come strong enough to accelerate protons from the dense 
parts of the thin shell to the dilute ones. 

The impact of this secondary ambipolar electric field 
on the protons in the thin shell is evidenced by Fig. 


FIG. 5: (Golor online) The thin shell at the time t = 68: 
Panel (a) shows np{x,y). Panel (b) shows {E^ -h Ey)^^‘^. 

which corresponds to the simulation time t = 6S. Figure 
ia) reveals a complicated proton density distribution 
within the thin shell. Most protons are still located be¬ 
hind the concave shocks at {x^y) = (20,250), (x^y) = 
(115,250), {x,y) = (70,350) and {x,y) = (160,350). 
New proton density maxima have appeared, which are 
centered at the positions x ^ 300 and y = 0, 45, 90 and 
135. These ^-coordinates correspond to those locations 
where the initial contact boundary xsiy) = 0 and where 
the density of the thin shell was lowest in Fig. [^a). The 
density decreases from 3 to about 1.8 near these maxima. 
These density jumps give rise to the strongest electric 
fields in Fig. [^b). The localized proton density deple¬ 
tions inside of the thin shell and their associated electric 
fields appear to be stable and they convect with the thin 
shell. This is evidenced by Movie 3 |33], which animates 

in time {E^ -h Ey)^^‘^. 

Figure demonstrates that two electrostatic shocks, 
which separate the downstream region at x ~ 300 from 
the upstream regions at large and low values of x, are 
forming at the time t = 68. The bulk of the protons is still 
propagating with the interpenetrating ion clouds with 
Vx ^ 0 or Vx ^ vq. a significant population is, however, 
also filling up the phase space interval with x ^ 300 and 
'^x ~ vq/2. The ion phase space structures at —100 < 
X < 200 and Vx ^ 0 and 400 < x < 500 and Vx ~ vq are 
fed by the shock-reflected ions and by ions that crossed 
the downstream region and were re-accelerated by the 
second shock. 

Summary: We have tested with a two-dimensional 
particle-in-cell simulation the stability of a thin shell of 
collisionless plasma, which is enclosed by plasma shocks. 
The thin shell was displaced along the normal of its av¬ 
erage plane. The displacement varied sinusoidally as a 
function of the position along the thin thell. We have 
demonstrated that the thin shell is unstable. The am- 
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FIG. 6: (Color online) The proton phase space density dis¬ 
tribution fi{x,y,px) at the time t = 68: The color range is 
linear. The x-position is given in units of IOOAd and Vx is 
given in units of vq. The full width of the y-axis is displayed. 


plitude of the displacement grew 15-fold. The thin shell 


followed a piecewise linear spatial distribution when the 
thin-shell instability saturated. The continuing accumu¬ 
lation of protons at the density maxima led to the forma¬ 
tion of large-amplitude density modulations within the 
thin shell. These proton density depletions, which are 
sustained by strong electric fields, will probably become 
phase space structures known as proton phase space holes 
m i3o]. The deflection of the upstream electrons and 
protons, which crossed the plasma shock, into opposite 
directions resulted in a net current within the thin shell 
that generated a magnetic field. 
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